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ABSTRACT 

Carbon monoxide (CO) is widely used as a tracer of molecular hydrogen (H2) in metal- 
rich galaxies, but is known to become ineffective in low metallicity dwarf galaxies. 
Atomic carbon has been suggested as a superior tracer of H2 in these metal-poor 
systems, but its suitability remains unproven. To help us to assess how well atomic 
carbon traces H2 at low metallicity, we have performed a series of numerical simulations 
of turbulent molecular clouds that cover a wide range of different metallicities. Our 
simulations demonstrate that in star-forming clouds, the conversion factor between 
[Cl] emission and H2 mass, Aci, scales approximately as Aci oc Z“^. We recover a 
similar scaling for the CO-to-H2 conversion factor, AqO) but find that at this point in 
the evolution of the clouds, Aco is consistently smaller than Aci, by a factor of a few 
or more. We have also examined how Aci and Aco evolve with time. We find that 
Aci does not vary strongly with time, demonstrating that atomic carbon remains a 
good tracer of H2 in metal-poor systems even at times significantly before the onset of 
star formation. On the other hand, Aco varies very strongly with time in metal-poor 
clouds, showing that CO does not trace H2 well in starless clouds at low metallicity. 
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1 INTRODUCTION 

Giant molecular clouds (GMCs) play a central role in the life 
of our own and other galaxies, as they are the sites where 
essentially all star formation occurs. For this reason, the 
study of the formation and evolution of GMCs is a highly 
important topic. However, this study is hampered by the 
fact that the two main constituents of the gas making up 
a GMC - molecular hydrogen (H 2 ) and atomic helium - 
do not emit radiation at the low gas temperatures present 
in most GMCs. In order to study the molecular gas, it is 
therefore necessary to observe emission from less abundant 
chemical species that do emit at low temperatures and that 
are correlated with the H 2 distribution. 

One of the most popular such tracer species is carbon 
monoxide (CO). This has several advantages as a tracer of 
molecular gas. First, there is a lot of it in most GMCs: it 
is typically the second most abundant molecular species, 
after H 2 . Second, it emits readily at low temperatures: the 
energy separation between the rotational ground state {J ~ 
0) and the first rotationally excited level (J = 1) is only 
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5.5 K. In addition, the main chemical pathways that lead to 
the formation of CO in the interstellar medium all involve 
H 2 as a key ingredient at some point (see e.g. the detailed 
discussion of CO formation chemistry in Glover et al.|2010 |, 
and so we would not expect large amounts of CO to form in 
regions devoid of H 2 . 

Observationally, CO emission has indeed proved to be 
a good tracer of molecular gas within the Milky Way. A 
number of different studies, summarized in the recent review 


by Bolatto, Wolfire & Leroy (20131, have shown that on 


scales comparable to or larger than individual GMCs there is 
a good correlation between the velocity-integrated intensity 
in the J = 1 —>■ 0 line of ^^CO and the H 2 column density. 
The conversion factor from integrated intensity to column 
density - the so-called X-factor - is given approximately by 


Aco,MW = 2 X 10^° cm ^(Kkm/s) ^ 


( 1 ) 


for Milky Way GMCs, and similar values are also found 
for other nearby. Milky Way-like disk galaxies. Simulations 
of GMCs in environments similar to that of the local ISM 
typically yield values for Aco close to this observationally 


determined value (see e.g. 

Glover & Mac Low||2011 

Shetty 

let al.|2011a|bl|Narayanan et al.|20111|Narayanan & 

lopkins 
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2013). Its lack of significant variation from cloud to cloud 


can be understood as a consequence of the fact that the main 
cloud properties that determine Xco in metal-rich galaxies 
- the velocity dispersion, the temperature and the surface 
density of the gas - do not vary much from cloud to cloud 
(Narayanan fc Hopkins|[2013 ). 


However, in low metallicity systems such as dwarf galax¬ 
ies, CO is much less useful as a tracer of molecular gas. A 
growing number of observational studies have provided ev¬ 
idence that Aco increases substantially as the metallicity 


decreases (see e.g. 

Israel]] 1997 

Madden et al.|1997 

Leroy et 

aL|20071|20111|Gratier et al.|2010t|Bolatto et aL|20111|Smith 

et al.[[2012 

1 , owing to the reduction in the mean extinction 


merical simulations also suggest that the Aco may become 
strongly time-dependent at low metallicity UGlover fc C5 a| 


2012b see also Section |3.2| below). Furthermore, CO has 


proved to be extremely difficult to detect in very low metal¬ 
licity systems: at Z ^ 0.1 Zq, there are secure detections of 
CO in only two systems, the Wolf-Lundmark-Melotte and 
Sextans A dwarf galaxies (Elmegreen et al. 112013 Shi et al. 
[2015] ). 

In view of the difficulties involved with using CO as a 
tracer of molecular gas in metal-poor galaxies, it is natural to 
ask whether there are better tracers that we might use. One 
interesting possibility is neutral atomic carbon, C. This pro¬ 
duces two fine structure lines in its electronic ground state, 
the ^Pi —>■ ^Po transition at 609/rm and the ^P 2 —>■ ^Pi tran¬ 
sition at 370/im. Its potential as a tracer of molecular gas 


was championed by Papadopoulos, Thi & Viti (2004), who 


pointed out that although simple PDR models predict that 
it will be confined to a thin transition region at the edge 
of molecular clouds (see e.g. Langer|[l976 1, observations do 
not bear this out, and have long shown that C and CO emis¬ 


sion are surprisingly well correlated in local clouds | Frerking 
et al.||1989| [Little et al.||199^ [Schilke et al.|pi)95| [Kramer 


et al. 2008). This behaviour can be understood as a conse¬ 


quence of clumpy substructure within molecular clouds: CO 
emission comes primarily from within clumps and filaments, 
while [C l] emission traces their surfaces, but since the typi¬ 
cal scale of the structures is small and there are often many 
along each line of sight, the result is that the projected CO 
and [C l] emission are well-correlated, both spatially and also 


in terms of velocity ( 

Stutzki et al.|1988| 

Genzel et al. 

1988 

Plume, Jaffe & Keene|1994[ [Kramer et al.|2004[ [Glover et al.[ 

2015 

). 


More recently, Offner et al. (2014) have explored the 


usefulness of atomic carbon as a molecular gas tracer in solar 


metallicity clouds. They used the 3d-pdr code (Bisbas et al. 


2012) to post-process simulations of turbulent molecular 


clouds, allowing them to predict the chemical composition 
and temperature of the gas at each point in the cloud. Based 
on this, they then computed synthetic emission maps for the 
609/im of [Cl] and the J = 1 —>■ 0 transition of CO, using 
the RADMC-3D radiation transfer code[^ They found that 
[C l] emission was indeed a good tracer of H 2 column den¬ 
sity in their model clouds. They also computed a mean [C l]- 
to-H 2 conversion factor, Aci, for their models, and showed 


^ http://www.ita.uni-heidelberg.de/~dullemond/software/radmc- 
3d/ 


that Aci = 1.1 x 10^^ cm ^(Kkm/s) ^ for a cloud with solar 
metallicity illuminated by the standard interstellar radiation 
field (ISRF), and Aci = 1.2 x 10^^ cm-2(Kkm/s)"^ for a 
solar metallicity cloud illuminated by a ten times stronger 


ISRF. However, Offner et al. (2014) did not explore whether 


Aci is sensitive to the metallicity of the gas. In addition, 
there are a couple of inconsistencies in the [Offner et al.[ 
(2014) approach that may influence their results. Firstly, 


they perform their hydrodynamical simulations assuming an 
isothermal equation of state ( Offner et al.[[20l3 l, and only 
solve for the actual temperature distribution during the later 
post-processing step. This is a reasonable approximation in 
dense, well-shielded gas, but not in the lower density, mod¬ 
erately shielded gas where atomic carbon dominates (c.f. 
their Figures 6 and 7). Secondly, they evolve the chemistry 
for a period of 10 Myr, but stop the hydrodynamical sim¬ 
ulation after ~ 1 Myr. They may therefore miss important 
non-equilibrium effects in the chemistry. 

In [Glover et al.[ (2015; hereafter. Paper I), we presented 
results from a study of the [C 1 ] emission produced by a 
model of a turbulent molecular cloud that does not suffer 
from these inconsistencies. In our simulation, we evolved the 
chemistry and internal energy of the gas simultaneously with 
the hydrodynamics, allowing us to capture any influence of 
the non-isothermal equation of state on the dynamical evo¬ 
lution and to ensure that the chemistry is not evolved for 
a longer period of time than the dynamics. The cloud that 
we modeled in Paper I had a solar metallicity, and was il¬ 


luminated by the standard Draine (1978) ultraviolet (UV) 


radiation field. We found, in agreement with [Offner et al.[ 
( [2014[ ), that the [Cl] emission line was a good tracer of the 
molecular mass of the cloud. 

In the present paper, we extend this study by examin¬ 
ing clouds with a range of different metallicities. For each 
metallicity, we run models with several different values for 
the strength of the ISRF. Our goal is to understand how 
well [C l] emission can be used as a tracer of molecular gas in 
this wider range of systems, and whether it becomes clearly 
superior to CO at low metallicities. Since clouds in extra- 
galactic systems are typically unresolved or only marginally 
resolved, our main focus in this study is on how Aci varies as 
a function of metallicity, radiation field strength and time. 

The structure of our paper is as follows. In Section 
we outline the numerical method used to perform our simu¬ 
lations and to post-process the results. We also give details 
of the initial conditions that we used for our simulations. 
In Section]^ we present and discuss our results. Finally, we 
conclude in Section (4] 


2 METHOD 

2.1 Hydrodynamical model 

To model the chemical, thermal and dynamical evolution of 
the gas in our model clouds, we use a modified version of the 


Gadget 2 SPH code (Springel 2005). Our modified version 


of the code includes a simplified chemical model that follows 
the abundances of H, H 2 , H+, C, C+, CO, HCO+ He, He+, 
and free electrons (Glover & Mac Low 2007a|b Nelson & 
Langer[[l999 ), a detailed cooling function, and a treatment 
of the effects of attenuation of the interstellar radiation field 


© 0000 RAS, MNRAS 000, 000-000 



















































































Is atomic carbon a good tracer of molecular gas in metal-poor galaxies? 3 


(ISRF) based on the treecol algorithm (Clark, Glover & 


Klessen 20121. To model gravitationally collapsing regions 


that become too small to resolve, we use a sink particle al¬ 


gorithm based on Bate, Bonnell & Price (19951 and Jappsen 


et al. (20051. The version of the code used here is the same 


as that used in Paper I, and a more extensive discussion of 
its capabilities can be found in that paper. 

Although our use of a simplified chemical network in¬ 
troduces some uncertainty into our results that would not 
be present if we used a more extensive network, the impact 
of this on our predicted C and CO abundances is likely to 


be small. We showed in an earlier paper (Glover & Clark 
|2012c| ) that the simplified chemical network that we use 
in this study produces comparable results for the C and 


CO distributions to the larger Glover et al. (20101 network. 


That network was in turn benchmarked against an exten¬ 
sive model derived from the UMIST astrochemical database 
( Woodall et al.|2007 1 and shown to produce results in good 
agreement with the latter for a range of different physical 
conditions. In addition, the fact that our results agree well 


with those of Offner et al. (20141, who make use of a cheml 


cal network containing over 200 species and several thousand 
reactions, also gives us confidence that our use of a simplified 
network is not a major source of error. 


2.2 Post-processing 

We construct synthetic emission maps of our simulated 
clouds in the —>■ ®Po fine structure line of [C l] (here¬ 
after simply referred to as the [C l] 1 —>■ 0 line) and the 
J = 1 —>■ 0 line of '^^CO using the radmc-3d radiation 
transfer code. The region modelled with radmc-3d is a cu¬ 
bical sub-volume of the cloud, with a side length of 16.2 pc 
that contains almost all of the CO and neutral atomic car¬ 
bon present in the model clouds. Before feeding our data to 
RADMC-3D, we first interpolate it onto a Cartesian grid with 
a resolution of 256^ zones. In Paper I we explored the effect 
of varying this resolution and showed that 256® zones was 
sufficient to fully resolve the [C l] emission. 

To compute the level populations required by the code, 
we used the large velocity gradient (LVG) approximation. 
The implementation of this approximation within radmc- 
3d is described in detail in Shetty et al. ( 2011a|b' I. We fixed 


the velocity resolution of the code at Av = 0.094km s“^ and 
accounted for small-scale unresolved velocity fluctuations 
by including a microturbulent contribution to the linewidth 
with size Umtrb = 0.2 kms“®. However, the values that we 
derived for the [Cl] and CO velocity-integrated intensities 
are relatively insensitive to these choices. 


2.3 Initial conditions 

Our default initial conditions in this study are similar to 
those u sed in paper I and in sev eral other previous papers 
2012a bl. The cloud mass is lO'^M©, 


(see e.g. Glover & Clark 


and the cloud is simulated using two million SPH particles. 
The resulting mass resolution is 0.5 M©. In paper I, we ex¬ 
plored the effects of increasing this resolution, but found 
that this resulted in only minor differences in the [C l] emis¬ 
sion produced by the cloud. The cloud was initially spherical, 
with a radius of 6.3 pc. The initial H nuclei number density 


of the gas was no = 276 cm“®, within a factor of two of the 
median value determined for local GMCs by |Roman-DuvaT| 
et al. (20101. The initial gas temperature was Tgas = 20 K 
and the initial dust temperature was Tdust = 10 kQ 


The velocity field of the gas in the clouds was con¬ 
structed as a Gaussian random field, with a power spectrum 
P{k) oc where k is the wavenumber. The amplitude of 
the velocity fluctuations was constrained to yield an initial 
RMS velocity of Urms = 2.8 kms“®. This value was chosen 
so that the initial kinetic energy of the cloud was equal to 
its gravitational binding energy. The velocity field was then 
allowed to evolve during the simulation without any further 
kinetic energy input. The cloud was simulated using vac¬ 
uum boundary conditions, but the artificial “evaporation” 
of material from the edge of the cloud was prevented by the 
application of an external pressure term (Benz|19OT|. 


We performed a series of simulations using these ini¬ 
tial conditions, and varying the metallicity of the gas and 
the strength of the imposed radiation field. The values used 
in the different simulations are summarized in Table |T] In 
our solar metallicity simulations, the total abundances of 
carbon and oxygen relative to hydrogen were taken to be 


xc = 1-4 X 10 and xq = 3.2 x 10 , respectively (Sem- 


bach et al.|2000 (. In the lower metallicity simulations, these 
values were scaled by a factor Z/Zq. In most of our simula¬ 
tions, we also assumed that the dust-to-gas ratio TO scaled 
as Z/Zq. However, there is some observational evidence that 
it actually falls off more rapidly as the metallicity decreases 


(e.g. Galametz et al. 2011 Herrera-Camus et al. 2012[ Remy 


Ruyer et al. 20141, and so we also performed three simu¬ 


lations (Z05-G1-D0375, Z02-G1-D01 and ZOl-Gl-DOOl) in 
which we adopted smaller dust-to-gas ratios. For simplic¬ 
ity, we neglect any effects arising due to changes in the size 
distribution or composition of the dust as we change the 
metallicity and dust-to-gas ratio. In all of our simulations, 
the oxygen and hydrogen were assumed to begin in neutral 
atomic form, and the carbon was assumed to start as C^. 


The spectral shape of our adopted interstellar radiation 


field follows Draine (19781 in the ultraviolet and Mathis, 
[Mezger fc Panagia| ( |198^ at longer wavelengths. We con¬ 
sider several different radiation field strengths, quantified by 
the parameter Go which is normalized so that Go = 1 corre¬ 
sponds to the |Draine| ( |1978] ) radiation field. In our treecol- 
based treatment of the attenuation of the radiation field, we 
used a total of 48 pixels to represent the angular variation 
of the radiation field around each SPH particle. In paper I, 
we verified that this was sufficient to yield well-converged re¬ 
sults. In most of our runs, we adopted a cosmic ray ionization 
rate for atomic hydrogen given by ((h = Ch.o = 10“®^ s“®. 
However, we also carried out a few simulations with larger 
values of ("h. Cosmic ray ionization rates for H 2 , CO etc. were 
computed by scaling the H ionization rates, as described in 


Glover & Clark (2012bI. 
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Table 1. List of simulations Table 2. Onset of star formation in each simulation 


ID 

Z/Z0 

Go 

Ch/Ch,o 

Notes 

ID 

tSF (Myr) 

Zl-Gl 

1.0 

1 

1 


Zl-Gl 

1.98 

Zl-GlO 

1.0 

10 

1 


Zl-GlO 

2.17 

Zl-GlOO 

1.0 

100 

1 


Zl-GlOO 

2.05 

Z05-G1 

0.5 

1 

1 


Z05-G1 

2.43 

Z05-G10 

0.5 

10 

1 


Z05-G10 

2.58 

Z05-G100 

0.5 

100 

1 


Z05-G100 

2.79 

Z03-G1 

0.3 

1 

1 


Z03-G1 

2.72 

Z03-G10 

0.3 

10 

1 


Z03-G10 

2.86 

Z03-G100 

0.3 

100 

1 


Z03-G100 

3.24 

Z02-G1 

0.2 

1 

1 


Z02-G1 

2.81 

Z02-G10 

0.2 

10 

1 


Z02-G10 

2.97 

Z02-G100 

0.2 

100 

1 


Z02-G100 

3.59 

ZOl-Gl 

0.1 

1 

1 


ZOl-Gl 

2.95 

ZOI-GIO 

0.1 

10 

1 


ZOI-GIO 

3.22 

ZOl-GlOO 

0.1 

100 

1 


ZOl-GlOO 

4.25 

Z003-G1 

0.03 

1 

1 


Z003-G1 

3.20 

Z003-G10 

0.03 

10 

1 


Z003-G10 

3.61 

Z003-G100 

0.03 

100 

1 


Z003-G100 

4.81 

Z05-G1-D0375 

0.5 

1 

1 

V = 0.375 X>0 

Z05-G1-D0375 

2.54 

Z02-G1-D01 

0.2 

1 

1 

X> = 0.1220 

Z02-G1-D01 

2.93 

ZOl-Gl-DOOl 

0.1 

1 

1 

22 = 0.01220 

ZOl-Gl-DOOl 

2.97 

ZI-GIO-CRIO 

1.0 

10 

10 


Zl-GlO-CRlO 

2.17 

ZI-GIOO-CRIOO 

1.0 

100 

100 


Zl-GlOO-CRlOO 

2.03 

Z03-G10-GR10 

0.3 

10 

10 


Z03-G10-GR10 

2.88 

Z03-G100-GR100 

0.3 

100 

100 


Z03-G100-CR100 

2.58 

ZOI-GIO-CRIO 

0.1 

10 

10 


ZOl-GlO-CRlO 

3.25 

ZOI-GIOO-GRIOO 

0.1 

100 

100 


ZOI-GIOO-CRIOO 

— 


D denotes the dust-to-gas ratio, with Dq being the value in the 
local ISM. Go is the strength of the ISRF relative to |Drain^ | |1978[ l 
and [Mathis, Mezger fc Panagia] ([1983^ . 


3 RESULTS 

3.1 The [C l]-to-H 2 conversion factor at the onset 
of star formation 

In Paper I we showed that in a typical turbulent Galac¬ 
tic cloud, the integrated intensity of the [Cl] 1 —>■ 0 line 
(VPci.i-o) scales linearly with the H 2 column density {N 112 ) 
only over a limited range of H 2 column densities. At high 
column densities, the correlation breaks down because the 
1 —>■ 0 line becomes optically thick, while at low column 
densities, the breakdown reflects the fact that most of the 
carbon in regions of low extinction is photoionized to by 
the external radiation field. Observations of [C l] in Galactic 
molecular clouds recover similar behaviour (see e.g. Beuther 
et al.lpoTd I. In the present paper, we show that increasing 
the radiation field strength and/or decreasing the metallic- 
ity does not substantially change this basic conclusion. This 
is illustrated in Figure]^ where we show two examples of the 
relationship between ITci,i-o and Nh 2 , taken from runs Zl- 
GIO and Z02-G1. In both cases, the values shown are for a 


time immediately before the onset of star formation, t = tsF- 
The actual physical time that this corresponds to for each 
cloud is listed in Table We see that in both cases, there 
is an approximately linear correlation between Wci.i-o and 
N 112 only over a relatively narrow range in H 2 column densi¬ 
ties. At low column densities, Wci,i-o falls off rapidly with 
decreasing AfH 2 , reflecting the fact that most of the carbon 
is photoionized to G”*" in regions with low dust extinction. 
At high column densities, on the other hand, the correlation 
between Wci,i-o and A^H 2 becomes sub-linear because the 
[Cl] 1 —>■ 0 line becomes optically thick. Moreover, even in 
the range of H 2 column densities where the mean value of 
IFci,i-o does scale approximately linearly with A’h 2 , we see 
that there is considerable scatter. 

It is clear from this that on small scales, there is no 
unique [C l]-to-H 2 conversion factor. Instead, we would ex¬ 
pect the conversion factor to vary significantly with position 
within the cloud, similar to the behaviour that we find for 


Ago (see e.g. Pineda, Caselli & Goodman 2008 Shetty et al. 


2011a|b Lee et al.|2014 1. A more interesting quantity to ex¬ 

amine is therefore the mean [C l]-to-H 2 conversion factor for 
the entire cloud. We define this as 


Y _ (AH2) 

Aci = 


(Wcia^o) ’ 


( 2 ) 


^ Note that since Tgas and Tiust rapidly adjust to close to their 
thermal equilibrium values following the beginning of the simula¬ 
tions, our results do not depend on these starting values. 


where the angle brackets denote the fact that we are averag¬ 
ing over the projected spatial extent of the cloud. We define 
Aci for the cloud in this way, rather than by computing its 
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[CI] 1-0 



column density [cm 


10 ^“ 10 ^“ 0.001 
Fraction of cloud area 



. 1 l: _I_ r ■i-imm Ktx-- ._,. . . 

1020 jq22 


Hg column density [cm 


10 "® 10 "® 10 "'* 0.001 
Fraction of cloud area 

Figure 1. (a) Two-dimensional PDF of the integrated intensity of 
the [C l] 1 —i 0 line, plotted as a function of the H 2 column density, 
, for run Zl-GlO at t = igp. The dashed line indicates a linear 
relationship between lFci,l-0 ^H 2 included to help 

guide the eye. The solid line shows the geometric mean of FFci,i—0 
as a function of AfH 2 • The points are colour-coded according to 
the fraction of the total cloud area that they represent, (b) As 
(a), but for run Z02-G1 at t = igp- 


value for each individual pixel and then averaging, for two 
reasons. First, this better reflects what we are actually doing 
when we make observations of unresolved molecular clouds 
and hence measure only the total emission rather than its 
spatial distribution. Second, this procedure allows us to deal 
sensibly with lines-of-sight that have H 2 without any signif¬ 
icant associated [C l] emission. 

When computing (A^Hs) and {Wci,i-»o), we include all 
of the pixels in our projections that have non-zero H 2 column 
densities. However, we note that this is not the only possi¬ 
ble way in which to compute ‘average’ H 2 column densities 
and [C l] integrated intensities for the simulated clouds. A 
commonly-used alternative method involves defining some 
minimum column density or intensity threshold, and then 
only averaging over the set of pixels with values above this 
threshold (see e.g. Clark & Glover 2015| or Sziics et ah. 


in prep.) This procedure better matches what is done in re¬ 
solved observations of clouds (e.g. Lee et al.|2014 1, while the 
procedure that we adopt is more appropriate if the clouds 
are unresolved. 

In Figure we show how Xci varies as a function 
of metallicity and radiation field strength for our simulated 
clouds. Once again, we show the results for t = tsF- The 
values of Xci shown in this hgure are given in terms of the 
canonical Galactic CO-to-H 2 conversion factor, Aqo.mw- 

Two points are immediately obvious. First, at so¬ 
lar metallicity, Xci is signihcantly larger than Xco.mw- 
In run Zl-Gl, we recover a value of Xci = 8.8 x 
10^° cm“^ km~^ s, within 25% of the value computed 
by Offner et al. (20141 for a somewhat denser cloud. This is 
more than a factor of 4 larger than Aco.mw, and roughly 
three times larger than the value of Xco that we recover 
for the same simulated cloud (see Figure]^). Increasing the 
strength of the ISRF by a factor of 10-100 increases Xci by 
around 30-50%, but also increases Aco by roughly the same 
amount. Therefore, at solar metallicity and at this point in 
the clouds’ evolution, there does not appear to be a com¬ 
pelling reason to prefer [C l] over GO as a tracer of the total 
cloud mass. 

Second, Figure demonstrates that Aci systemati¬ 
cally increases as we reduce the metallicity of the gas, with a 
scaling that is close to Aci oc in the Go = 1 runs. In the 
runs with stronger radiation fields, we see a slightly steeper 
scaling for metallicities in the range 0.2 < Z/Zq < 1.0, 
which then flattens as we move to even lower metallicities. 

For comparison, we illustrate in Figure how the 
value of Ago that we measure at the onset of star forma¬ 
tion varies with metallicity. As in the case of Aqi, we see 
that at Z is 0.2 Z©, Ago increases as the metallicity de¬ 
creases, with a dependence that again scales roughly as Z~^. 
At lower metallicities, the behaviour becomes more com¬ 
plicated. When Go = 1, Ago continues to increase with 
decreasing metallicity, but in the Go = 10 and Go = 100 
runs we find instead that Ago remains constant or even 
decreases with decreasing Z. This occurs because in these 
conditions, the total CO emission of the cloud becomes 
dominated by the contribution from a single dense, com¬ 
pact, gravitationally-collapsing clump (see Figure [^. The 
high column density of the gas in this clump provides ef¬ 
fective shielding for the CO even when the metallicity is 
very low, and allows the gas to build up a large CO column 
density within the clump. Consequently, the CO emission 
from the clump becomes optically thick, at which point it 
becomes only weakly sensitive to the fractional abundance 
of CO within the dense gas. As a result, the mean value of 
the ^^CO integrated intensity, (VFgo), which varies strongly 
with metallicity at high Z, starts to vary only weakly with 
metallicity once we reach this clump-dominated regime. At 
the same time, however, the H 2 content of the cloud still de¬ 
creases substantially with decreasing metallicity even once 
we reach this clump-dominated regime, since even at the 
lowest metallicities, much of the H 2 is still located outside 
of the dense clump. 

The results of our comparison show that at this point 
in the evolution of the cloud, the cloud-averaged CO emis¬ 
sion remains a better tracer of the molecular mass than the 
cloud-averaged [C l] emission, even at metallicities as low as 
Z = 0.03 Zq. This can be seen more easily if we plot the 
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Figure 2. (a) [Cl]-to-H 2 conversion factor at t = igp, computed 
for the [C l] 609 lira line, and plotted in units of the Galactic CO- 
to-H 2 conversion factor mw = 2 x cm“^(Kkm/s)“^) 

for a range of different metallicities and radiation field strengths. 
We show results for radiation field strengths Go = 1 (solid line), 
Go = 10 (dotted line) and Go = 100 (dash-dotted line). The 
dashed line indicates a slope of and is included merely to 
guide the eye. (b) As (a), but for the CO-to-H 2 conversion factor, 
-’i'co- 


ratio of Xco to Xci, which is the same as the ratio of the 
cloud-averaged [C l] and CO integrated intensities. This is 
shown in Figure We see that in clouds at the onset of 
star formation, {Wci)/{Wco) decreases as we decrease the 
metallicity, although the dependence is fairly weak: an order 
of magnitude drop in Z leads to a decrease in the intensity 
ratio of no more than a factor of two. We also see that except 
at the very lowest metallicity that we study, {Wci)/{Wco) 
is highly insensitive to the value of Gq. 


3.2 Time evolution of Xci and Xqo 

In the previous section, we examined the behaviour of Xci 
and Xco in the case where the gas has already collapsed 



X [pc] 


Figure 3. (a) Map of the integrated intensity of the J = 1 —> 0 
transition of in run Zl-Gl at time t = tgp. (b) As (a), but 

for run ZOI-GIOO. 


and started forming stars, and showed that in this case, CO 
is a superior tracer of the molecular mass, at least for the 
range of metallicities and UV field strengths considered here. 
However, it does not immediately follow that CO will be a 
superior tracer of H 2 at all times during the evolution of the 
cloud. As we have already seen, at low metallicity the bulk of 
the CO emission produced by the cloud comes from between 
one and a few dense, compact, collapsing clumps (see also 


Glover & Clark 2012b who report similar results). These 


clumps have higher volume densities and higher extinctions 
than the bulk of the cloud, and hence are the regions in 
which CO is best able to resist photodissociation. However, 
the dynamical timescale on which these clumps evolve is rel¬ 
atively short: they will either collapse on a free-fall timescale 
{tg ~ 0.5 Myr for clumps with a mean density n ~ 10"^cm“^, 
or even shorter at higher densities), or will be disrupted by 
turbulence on a similar timescale. This prompts us to ask 
how Xci and Xco evolve as a function of time: do we find 
larger values at earlier times, when these dense clumps have 
not yet formed? 
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Figure 4. Ratio of the mean integrated intensities (in K km s“^) 
of the [Cl] 609 fim line ((Wci)) and the J = 1 —>■ 0 line 

((kFco})i computed at t = tgp and plotted as a function of metal- i— 

licity. We show results for radiation field strengths Go = 1 (solid 1—. 

line), Go = 10 (dotted line) and Go = 100 (dash-dotted line). 

6 

To investigate this, we have calculated Xci and Xco „ 

for each of our clouds at a time t = tsF — 1 Myr. The results g 

are plotted in Figurej^ We see that when Go = 1, the values “ 

of Xci that we measure at t = tsF — 1 Myr are similar to 'o 

or smaller than those that we measure at t = tsF- In other x 

cv 

words, in these conditions, [C l] emission remains a good ^ 

tracer of the molecular mass of the cloud even prior to the >< 

onset of gravitational collapse. The same is true in the runs 
with Go = 10 and 100 at solar metallicity, but not at lower 
metallicities, where we find a much steeper increase in Xci 
with decreasing Z than we did in the previous section. That 
said, the actual increase in the values of Xci relative to those 
in the t = tsF case remains relatively modest: at most, we 
see increases of around a factor of a few. 

The story is very different in the case of CO. At solar 
metallicity, much of the CO emission comes from gas with a 
number density of a few thousand cm“^. Regions with this 
density are found throughout the cloud (see e.g. Figure 5 
in [Glover fc Clark||2012b| and most of this gas is not self- 
gravitating. Individual CO-bright regions are formed and de¬ 
stroyed by the action of the turbulence, but the total amount 
of gas in these regions does not vary much. Consequently, 

Xco does not vary strongly with time. At much lower metal¬ 
licities, however, dense, self-gravitating clumps dominate the 
total emission, and the value of Xco becomes highly sensi¬ 
tive to the evolutionary state of these clumps. In this case, 
the values that we measure for Xco at t = tsF — 1 Myr are 
much, much larger than those at t = tsF- 

This point is further emphasized if we look at the 
[Cl]/CO integrated intensity ratio at the same time (Fig¬ 
ure]^. At metallicities close to solar, the ratio of the [C ijand 
CO integrated intensities at this point in the cloud’s evolu¬ 
tion is around 0.5 to 0.6 (for values expressed in terms of 
K km s“^), somewhat larger than at t = tsF, but of the 
same order of magnitude. At lower metallicities, however, 
the ratio increases dramatically as the diffuse CO emission 



Z / Ze 


Figure 5. (a) As Figure]^, but showing values of Xci computed 
at t = tsF — 1 Myr. (b) As (a), but for Xco- 


vanishes. When Go = 1, this sharp increase in the [Cl]/CO 
intensity ratio occurs at a metallicity of around 0.1 Zq, but 
increasing Go leads to this change of behaviour happening 
at higher metallicity, owing to the more effective destruc¬ 
tion of the diffuse CO in runs with a higher radiation field 
strength. 

In Figure [^, we explore in more detail how Xci and 
Xco evolve with time in several representative runs. In run 
Zl-Gl, Xci rises gradually over time as the cloud collapses 
and an increasing fraction of the available carbon is con¬ 
verted into CO. In run Z02-G10, on the other hand, the 
collapse of the cloud leads to a decrease in Xci, as in this 
case, the increased conversion of C"*" to C is more important 
than the conversion of C to CO. Finally, in run Z003-G100, 
there is no clear trend in the evolution of Xci with time. 
Importantly, we see that in no case does Xci vary by more 
than a factor of a few over the time period examined, which 
corresponds to roughly half a free-fall time. 

If we now look at the evolution of Xco in the same 
three runs over the same period of time (Figure Iz^) , we 
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Figure 6. As Figure]^ but for t = fgF ~ 1 Myr. We show results 
for radiation field strengths Go = 1 (solid line), Go = 10 (dotted 
line) and Go = 100 (dash-dotted line). Note the difference in 
vertical scale between this figure and Figure]^ 

see that the behaviour is very different. In run Zl-Gl, the 
dense clumps and filaments created by the turbulence pro¬ 
vide enough shielding to allow significant quantities of CO 
to form even prior to the onset of run-away gravitational col¬ 
lapse. Consequently, the CO content of the cloud does not 
change substantially over time and Xco remains approxi¬ 
mately constant. 

In run Z02-G10, on the other hand, the combination of 
the lower dust extinction and higher radiation field strength 
mean that the over-dense structures created by the turbu¬ 
lence alone contain little CO. In this case, it is the gravi¬ 
tational collapse of some of these structures that triggers a 
substantial increase in the CO luminosity of the cloud. The 
result is that in this run, Xco evolves strongly over time, 
decreasing by almost a factor of 100 over the 1.5 Myr period 
shown here. 

Finally, we see that in run Z003-G100, the time- 
dependence of Xco becomes even more extreme. In this 
run, CO forms only at very high densities, shortly before the 
onset of star formation, and so Xco rapidly changes from 
being orders of magnitude larger than the Galactic value at 
t tsF to only 5-10 times larger than the Galactic value at 
t ~ tsF and later times. 

3.3 Dependence of Wci and Wco on extinction 

To help us understand why Xci is much less time dependent 
than Xco in low metallicity clouds, it is useful to look at 
how the integrated intensities of the [C l] 609 /rm line and 
the CO J = 1 —>■ 0 line depend on the line-of-sight visual 
extinction of the gas and how this dependence changes as 
we vary the metallicity. In Figure [^, we show a cumulative 
plot of the fraction of the total [C l] emission coming from 
lines-of-sight at or below the specified visual extinction for 
a number of different simulations that cover the full range 
of metallicities that we examine in this paper. In each case, 
the curves correspond to the state of the cloud at f = tsp- 



t - [Myr] 


Figure 7. (a) Evolution with time of Xci in runs Zl-Gl (solid 
line), Z02-G10 (dashed line) and Z003-G100 (dotted line), (b) As 
(a), but showing the evolution with time of Xco- 


In Figure [^, we show a similar plot for the CO J = 1 —^ 0 
emission in the same set of simulations. 

We see in the solar metallicity run, almost all of the 
[Cl] emission is confined to lines-of-sight with 1 < Av < 10, 
with roughly 50% coming from the narrower range 2.5 < 
Av < 5. However, as we decrease the metallicity, the range 
of extinctions occupied by the lines-of-sight producing most 
of the [C l] emission shifts to systematically lower values, so 
that by the time we reach the 0.03 Zq run, almost all of the 
emission is coming from gas with Av < 1. The amount by 
which the curve shifts from run to run roughly corresponds 
to the change in the dust-to-gas ratio, which is consistent 
with the same range of gas column densities dominating the 
emission in each case, regardless of the dust-to-gas ratio. 
This range of column densities extends from a few times 
lO^'^ cm“^ to ~ 10^^ cm“^, and it is important to note that 
the mean H nuclei column density of the cloud at t = 0, 
which is approximately (XH)init — 7 x 10^^ cm~^, lies in the 
middle of this range. It is therefore not surprising that the 
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Visual extinction 


Figure 8. (a) Fraction of the total [C l] emission from the sim¬ 
ulated clouds coming from lines of sight with values of Ay at or 
below the specified value, measured at f = tgp. Results are plot¬ 
ted for runs Zl-Gl (solid), Z05-G1 (dotted), Z03-G1 (dashed), 
Z02-G1 (dot-dashed), ZOl-Gl (dot-dot-dot-dashed) and Z003-G1 
(long-dashed), (b) As (a), but for the total J = 1 —>■ 0 

emission. 


clouds are bright in [C l] right from the start of the simula¬ 
tion, and hence that Xci does not vary strongly with time. 

In the case of CO, however, we see quite different be¬ 
haviour, as Figure demonstrates. In run Zl-Gl, most 
of the CO emission is produced by lines of sight with 
2 < Ay < 5, similar to the behaviour of the [C l] emis¬ 
sion. As we decrease the metallicity, we again find that the 
range of extinctions responsible for producing most of the 
emission shifts to lower values. However, in this case, the 
lines shift by a much smaller amount: a factor of 30 change 
in metallicity and dust-to-gas ratio shifts the curves by no 
more than around a factor of three. This means that the 
range of gas column densities that dominate the total CO 
emission shifts to higher values as we decrease the metallic¬ 
ity. For example, in run ZOl-Gl, around 50% of the total 


CO J = 1 —f 0 emission is produced by lines-of-sight with 
Ay > 1, corresponding to A^h — 2 x 10^^ cm~^. 

This behaviour provides a simple explanation for the 
strong time dependence that we see in Xco at low metal- 
licities. In our low metallicity clouds, the column density 
required to provide sufficient dust shielding to allow CO to 
form in significant amounts is much higher than the mean 
cloud column density. We know that in clouds dominated 
by turbulence, the column density PDF is log-normal, with 
a typical dispersion around the mean column density of no 
more than a factor of two in clouds in which the turbulence 


is primarily solenoidal in nature (see e.g. Federrath et al. 


2010 Burkhart fc Lazarian 201^[Abreu- Vicente et al. 20151. 
Consequently, in turbulence-dominated clouds prior to the 
onset of gravitational collapse, there are few or no lines-of- 
sight associated with bright CO emission. It follows from 
this that the mean CO brightness of the cloud is very small, 
and hence that Xco is very large, since many of these clouds 
still contain significant amounts of H 2 . However, once run¬ 
away gravitational collapse begins, the column density PDF 
typically develops a pronounced power-law tail at high col¬ 


umn densities (see e.g. Klessen 2000 Kainulainen et al. 2009 


Schneider et al. 20151. This increases the number of lines- 


of-sight that have high extinctions and hence significant CO 
brightnesses. The result is a dramatic increase in the mean 
CO brightness of the cloud, which drives a dramatic de¬ 
crease in Xco- Note that we do not see this behaviour in 
high metallicity clouds because in that case, the mean cloud 
column density already lies in the range of column densi¬ 
ties for which we find bright CO emission. Therefore, in this 
case, the cloud is CO-bright prior to the onset of runaway 
gravitational collapse, and so Xco does not vary strongly 
with time. 


3.4 Effect of varying the dust-to-gas ratio 

In most of our simulations, we assume that the dust-to¬ 
gas ratio TO is linearly proportional to the metallicity Z/Zq. 
Observationally, this seems to be a reasonable assumption 
for metal-rich galaxies like our own Milky Way ( |Sandstrom| 
et al. 20131. However, there is growing evidence that in 


metal-poor systems, such as dwarf galaxies, the dust-to¬ 
gas ratio becomes a much steeper function of the gas-phase 


metallicity (see e.g. Galametz et al. 2011 Herrera-Camus 
et al. 2012 Remy-Ruyer et al.||2014l. It is therefore useful 


to examine what effect this steeper dependence of on Z 
may have on our results for Xci and Aco- 

We have therefore performed three additional simula¬ 
tions in which the dust-to-gas ratio was not simply scaled 
with Z/Zq. In run Z05-G1-D0375, we adopted a gas-phase 
metallicity of Z = 0.5 Zq and a dust-to-gas ratio oi TO = 
0.3751?0. In run Z02-G1-D01, we adopted a gas-phase 
metallicity Z = 0.2 Zq and a dust-to-gas ratio TO = 0.1 Oq. 
Finally, in run ZOl-Gl-DOOl, we adopted a gas-phase metal¬ 
licity Z = 0.1 Zq and a dust-to-gas ratio TO — O.OIUq. The 
values in the first of these runs were chosen to be a good 
match for the gas-phase metallicity and dust-to-gas ratio in 
the Large Magellanic Cloud ( Gordon et al.|2014 l, while the 
values in the other two runs were chosen to be consistent 
with the relationship between metallicity and dust-to-gas 
ratio inferred by |Remy-Ruyer et al. (20141. 

In Table we list the values of Aci and Xco in these 
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Table 3. Comparison of Xci and Xqo in runs with low and 
standard dust-to-gas ratios 


ID 

Xci(ti) 

Xco(ti) 

Xci(tSF) 

Xco(tSF) 

Z05-G1-D0375 

4.21 

2.21 

7.15 

1.76 

Z05-G1 

3.85 

1.79 

7.05 

1.69 

Z02-G1-D01 

18.7 

35.1 

16.2 

3.57 

Z02-G1 

18.2 

10.9 

18.6 

3.13 

ZOl-Gl-DOOl 

9.32 

3420 

10.5 

6.15 

ZOl-Gl 

54.0 

183 

40.9 

7.75 


Xqi and Xqo nre given in units of the Galactic CO-to-H 2 con¬ 
version factor i^co.MW (s®® E*!' 0- Values are shown for tgF 
and = tsF ~ 1 Myr. 


runs at tsp and ti = tsF — 1 Myr. For comparison, we also 
list the values of Xci and Xco at the same times in runs 
with the same metallicity and radiation field strength but 
with a dust-to-gas ratio scaled linearly with metallicity. 

We see from the table that in run Z05-G1-D0375, the 
lower dust-to-gas ratio has only a small effect on the outcome 
of the simulations. We recover systematically higher values 
of Xci and Xco than in run Z05-G1, but only by a small 
amount: around 10-20% at ti, and less than 5% at f = fsF- 

In run Z02-G1-D01, we also find only minor differences 
compared to run Z02-G1 if we look at Xci or at the value 
of Xco measured at tsF- However, we see that there is a 
much larger difference in Xco at earlier times: at t = fi, it 
is roughly three times larger in run Z02-G1-D01 than in run 
Z02-G1. This difference in behaviour is again a consequence 
of the fact that at f = fi, the CO emission from this cloud 
is dominated by emission from diffuse gas, while at t = tsF, 
it is dominated by emission from a single dense core. The 
visual extinction of the dense core is sufficiently high that 
even if we decrease the dust-to-gas ratio by a factor of two, it 
still remains well-shielded, allowing CO to resist photodis¬ 
sociation. Consequently, the change in V does not have a 
strong influence on the value of Xco that we recover at this 
point in the evolution of the cloud. At earlier times, how¬ 
ever, most of the CO emission is produced by gas with a 
relatively low visual extinction. The amount of CO that can 
survive in this gas is much more sensitive to changes in Ay, 
and hence to change in the dust-to-gas ratio, and so it is not 
surprising that in this regime, Xco is highly sensitive to the 
value of T). 

Finally, in run ZOl-Gl-DOOl, we see that the substan¬ 
tial reduction in D that is inferred for very low metallicity 
systems leads to significant changes in both Xci and Xco- 
Reducing by a large amount - a factor of ten in this case - 
has two main effects. First, it significant reduces the amount 
of H 2 formed in the cloud, since the H 2 formation rate is 
proportional to "D. Second, it reduces the amount of dust 
shielding, which decreases the equilibrium abundances of 
both atomic carbon and CO. In the case of [Cl], the former 
effect dominates: although the mean [C 1 ] integrated inten¬ 
sity decreases, the H 2 column density decreases by a much 
larger factor, and so Xci decreases. For CO at early times 
(i.e. prior to runaway gravitational collapse), the reduction 
in shielding dominates instead. The CO abundance and CO 
integrated brightness both decrease by a much larger factor 
than the H 2 column density, and so Xco increases signifi¬ 


cantly. Finally, at t = tsF, when the total CO emission of 
the cloud is dominated by the contribution from the gravita¬ 
tionally collapsing core, Xco becomes largely independent 
of T), decreasing by only 20% for a factor of ten decrease in 

V. 

3.5 Effect of varying the cosmic ray ionization 
rate 

In most of the simulations presented in this paper, we have 
kept the cosmic ray ionization rate of atomic hydrogen fixed 
at its default value of ()h,o = s“^, even when vary¬ 

ing the strength of the interstellar radiation field. We have 
done this in order to allow us to focus solely on the effect 
of varying one parameter at a time. However, it is likely 
that this is not a good approximation for the behaviour 
of real star-forming systems. The low energy cosmic rays 
that make the dominant contribution to ("h are believed to 
be produced by Fermi acceleration in supernova remnants. 
Therefore, the cosmic ray energy density should scale lin¬ 
early with the supernova rate, provided that the mean time 
between supernovae is shorter than the timescale on which 
the cosmic rays escape from the galaxy. The supernova rate 
in turn depends on the star formation rate, averaged on a 
10-20 Myr timescale, and so it is reasonable to assume that 
the cosmic ray energy density, and hence the cosmic ray 
ionization rate, should scale linearly with the star formation 
rate (see e.g. |Papadopoulos|2010[[Papadopoulos et al.|201ip . 
The ultraviolet portion of the ISRF that is responsible for 
the photodissociation of H 2 and CO and the photoioniza¬ 
tion of atomic carbon is dominated by emission from mas¬ 
sive stars, and so should also scale linearly with the star 
formation rate. Therefore, in a real star-forming system, it 
is likely that increases in Go are associated with comparable 
increases in ((h- 

To investigate the effect that this has on the relative 
importance of [Cl] and CO as tracers of molecular gas, we 
have carried out several simulations in which we increased 
both Go and :Ch by factors of 10 and 100 compared to our 
default values. The details of these runs are listed at the 
end of Table In five of these runs, star formation be¬ 
gins after roughly 2-3 Myr, just as in the runs in which we 
keep (Ch fixed. In run ZOI-GIOO-CRIOO, however, the com¬ 
bination of the elevated ISRF and stronger cosmic ray flux 
provide so much heating that star formation is completely 
suppressed. The gas in the cloud heats up rapidly, with the 
mass-weighted mean temperature reaching T ~ 1000 K after 
only 1 Myr, and T ~ 2400 K after 2 Myr. This strong heat¬ 
ing causes the cloud to expand, reducing its density and 
allowing the temperature to climb further, a process that 
eventually results in the destruction of the cloud. As a re¬ 
sult of this, the C and CO abundances remain extremely 
small throughout the simulation. For this reason, we do not 
include this run in our further analysis. 

For the remainder of the runs with elevated cosmic ray 
ionization rates, we compute Xci and Xco at tsF and ti = 
tsF — 1 Myr. The results are shown in Table along with 
the values from the corresponding runs with fixed ^h. 

We see from Table that an increase in ("h leads to an 
increase in Xco and a decrease in Xci. This behaviour is 
easy to understand chemically. As we increase the cosmic 
ray ionization rate, we increase the amount of He^ present 


© 0000 RAS, MNRAS 000, 000-000 









Is atomic carbon a good tracer of molecular gas in metal-poor galaxies? 11 


Table 4. Comparison of Xqi 
cosmic ray ionization rates 

O 

O 

a 

in runs with 

different 

ID 

Xci(ii) 

Xco(tl) 

Xci(tSF) 

Xco(tSF) 

ZI-GIO-GRIO 

2.46 

2.04 

4.61 

2.01 

Zl-GlO 

2.95 

1.47 

5.85 

1.58 

ZI-GIOO-GRIOO 

2.87 

11.7 

3.18 

3.02 

Zl-GlOO 

4.94 

2.86 

6.85 

1.99 

Z03-G10-GR10 

36.3 

48.1 

20.9 

7.35 

Z03-G10 

54.5 

35.2 

33.7 

6.1 

Z03-G100-GR100 

85.3 

2110 

22.0 

20.8 

Z03-G100 

326 

810 

64.5 

15.7 

ZOI-GIO-GRIO 

213 

30500 

85.6 

25.1 

ZOI-GIO 

463 

34300 

132 

15.0 


Xqi and Xco given in units of the Galactic CO-to-H 2 con¬ 
version factor -^co,MW (s®® 0- Values are shown for tgp 

and t\ = tsF — 1 Myr. 


in the gas. This in turn increases the rate at which CO is 
destroyed by the dissociative charge transfer reaction 

CO + He+ ^ C+ + O + He, (3) 

which is the dominant destruction mechanism for CO in 
well-shielded gas. Increasing the cosmic ray ionization rate 
therefore reduces the amount of CO present in the cloud, 
which in turn reduces the mean CO brightness and increases 
Xco- In well-shielded regions, the C^ produced by this re¬ 
action does not persist in the gas for long. Instead, it is 
converted to neutral atomic carbon by radiative recombi¬ 
nation. The timescale for this process is much shorter than 
the timescale on which these carbon atoms are incorporated 
into new CO molecules, and so the end result is an increase 
in the C abundance, and hence an increase in Wci and de¬ 
crease in Xci- Note also that although our chemical model 
does not include the reactions 


C + He+ - 

and 

-4 C+ -f He, 

( 4 ) 

C + H+ - 

^ C+ + H, 

( 5 ) 


we do not expect these reactions to signihcantly affect the 
atomic carbon abundance, despite the increased H'*’ and 
He'*' abundances produced by the higher cosmic ray ioniza¬ 
tion rate. The reason for is that, in contrast to reaction 


which proceeds at close to the Langevin rate (Anicich & 
Huntress 19861, reactions and occur only very slowly, 


since both have rate coefficients that are roughly a million 
times smaller than that of reaction ( [Kimura et al.|[l993| 
Stancil et al. JT^. 


Looking in more detail at the values for Xci and Xco 
listed in Table we see that at t = tsF, even relatively 
large changes in the cosmic ray ionization rate lead to only 
comparatively small changes in the conversion factors. Xco 
increases by less than a factor of two, even if we increase ("h 
by a factor of 100. Xci is somewhat more sensitive to the 
value of but even in this case varies by no more than a 
factor of three for a factor of 100 change in the ionization 
rate. At earlier times in the evolution of the cloud, Xco 
becomes more sensitive to (^h and Xci becomes slightly less 



Figure 9. Ratio of the mean fractional abundances of CO and 
H 2 , plotted as a function of the density, for runs Zl-GlOO (solid 
line) and ZI-GIOO-CRIOO (dashed line) at t = tgp. Note that in 
fully molecular gas, our choice for the total carbon abundance 
implies that a;co/^H2 = 2.8 X 10“^. 


sensitive, but again the largest differences we see are around 
a factor of a few. Therefore, the qualitative details of the 
results presented in the previous sections are robust against 
changes in the cosmic ray ionization rate, even though the 
quantitative details are not. 

Finally, it is interesting to compare our results with the 
recent prediction by Bisbas, Papadopoulos, & Viti (20151 
that cosmic ray ionization rates of around 10“^°--10“^“ s“^ 
are sufficient to render molecular gas largely untraceable by 
CO at typical molecular cloud densities. As is clear from the 
results in Table|^ values of (^h in this range do increase Xco 
slightly in evolved clouds, but the clouds are still bright in 
CO, which remains a good tracer of molecular gas at metal- 
licities close to solar. Our results therefore do not conhrm 
the Bisbas, Papadopoulos, & Viti (20151 prediction. The dis¬ 
crepancy between our study and theirs seems to arise from 
the difference in how we treat the density structure of the 
model clouds. |Bisbas, Papadopoulos, fc Viti| largely focus 
on the effect of cosmic rays in uniform density clouds, and 
hnd that in this case, at typical GMCs densities of 100- 
1000 cm“®, CO is strongly suppressed when (^h (S lO”'^® s“^ 
even in the absence of photodissociation. In our turbulent 
clouds, we also find low CO abundances in gas with this 
range of densities when ("h is large, as illustrated in Figure]^ 
In this respect, our model produces results in good agree¬ 
ment with the [Bisbas, Papadopoulos, fc Viti] results. How¬ 
ever, thanks to the influence of the turbulent velocity field, 
the density distribution inside our model clouds is highly 
inhomogeneous, with a significant mass fraction located in 
regions with a mean density substantially higher than the 
volume-weighted mean density of the cloud. For example, 
in run Zl-GlOO-CRlOO at t = tsF, almost 10% of the mass 
of the cloud is located in clumps with n > 10^ cm“®, where 
both our model and the Bisbas, Papadopoulos, & Viti (20151 
model agree that the carbon will be almost entirely locked 
up in CO. These dense clumps therefore remain CO-bright. 
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Moreover, since these clumps are distributed throughout 
much of the cloud, the mean CO brightness of the cloud 
as a whole remains high, with the consequence that Xco 
remains close to its standard Galactic value. We can be con¬ 
fident that this is not simply due to some peculiarity of 
the model cloud studied here because our results for solar 
metallicity clouds are in reasonable agreement with those 
presented in the study of Clark & Glover (20151, which ex¬ 
amined a wider range of cloud masses and densities and used 
a different random realization of the turbulent velocity field. 
Instead, this seems to be a general consequence of the pres¬ 
ence of supersonic turbulence within the clouds. We note 
that the fact that the cloud inhomogeneity created by tur¬ 
bulence would allow more CO to survive was mentioned as a 
possible caveat by Bisbas, Papadopoulos, & Viti (20151, al¬ 
beit in the context of clouds with higher than normal levels 
of turbulence. Our results show that this effect is important 
even in clouds with a ‘normal’ level of turbulence (i.e. ones 
which sit on the standard size-linewidth relation, as derived 
by e.g. Larson||1981 or Solomon et al. 19871. 


3.6 Missing physics 

Our models neglect several physical processes that will be 
active in real molecular clouds. The most important of these 
are the freeze-out of CO onto cold dust grains, the influence 
of the magnetic field on the dynamical evolution of the gas, 
and the effects of star formation and the associated feedback 
on the chemical, thermal and dynamical state of the cloud. 
In this section, we briefly discuss how these missing processes 
may affect our results. 


3.6.1 Freeze-out 


cosmic-ray induced desorption due to local spot-heating of 
the grains. Our neglect of photodesorption, which is more 
difficult to model, again means that we will tend to overes¬ 
timate the extent to which CO is depleted. 

Having made these assumptions, we can relate the frac¬ 
tional abundance of gas-phase CO, a;co,gas to the fractional 
abundance of CO ice, a:co,ice using the following equation: 

I7acc3:cO,gas — (It(therm “t“ I?cr) 3:cO,ice. (6) 


Here, i?acc is the rate per CO molecule at which CO is ac¬ 
creted by dust grains, i?therm is the rate per CO molecule at 
which CO evaporates from the grain surfaces and Ra is the 
rate per CO molecule at which CO is desorbed by cosmic 
rays. Following Hollenbach et al. (20091, we have 

Ra. 




T>n s ^ molecule ^, 


where n is the number density of hydrogen nuclei. 


Rtherm = 1-04 X 10 exp 


m 


s ^ molecule 


( 7 ) 


( 8 ) 


( 9 ) 


where Td is the dust temperature, and 

Rcr = 6.0 X 10~^^Ch,i7 s~^ molecule”^, 

where Ch,i 7 = Ch/10~^^ s“^. The expression here for the 
thermal desorption rate assumes a CO adsorption energy of 
E-o/k = 960 K, following Aikawa et al. (19961. The cosmic 
ray desorption rate is based on the calculations of|Herbst fc| 


Cuppen (20061. 


Armed with the depletion and desorption rates, we can 
write the ratio of gas-phase CO to CO ice as 

^COjgas _ I^therm 4“ Rcr (10) 

XCO.ice R acc 

For Td ^ 17 K, cosmic-ray desorption dominates over ther¬ 
mal desorption and this reduces to 


In dense gas, the mean time between collisions between a 
given gas-phase particle (e.g. an H 2 or CO molecule) and a 
dust grain can become shorter than the dynamical time of 
the gas. In this regime, the freeze-out of gas phase chemi¬ 
cals onto the surface of dust grains can become important, 
provided that the average residence time of these species on 
the grain surface is also long. The residence time depends 
sensitively on the binding energy of the species to the sur¬ 
face, and so H 2 , which binds only very weakly, does not 
freeze-out, while CO, which binds much more strongly, can 
become significantly depleted from the gas phase (see e.g. 
Bergin et al.||20(i2 for a good observational example). 


The simplified chemical model that we use in our simu¬ 
lations does not account for the freeze-out of CO, and so the 
question arises as to whether this leads us to over-estimate 
the mean CO integrated intensity, and hence under-estimate 
Xco. To investigate this, we have attempted to account ap¬ 
proximately for the effects of CO freeze-out in an additional 
post-processing step. We make two important simplifying 
assumptions. First, we assume that locally, the ratio of gas- 
phase CO to CO ice has reached equilibrium. This is a rea¬ 
sonable assumption in dense gas, but means that we will 
overestimate the degree to which CO is depleted in diffuse 
gas, particularly at sub-solar metallicities. Second, we as¬ 
sume that the only processes that affect the ratio of gas- 
phase CO to CO ice are the accretion of CO molecules by 
dust grains, thermal desorption of the CO molecules and 


( 11 ) 


^CO,gas _ Rcr ....^ g ^ iQ 4 y— 1/2 Ch,17 

^CO,ice Race T)tI 

In dense gas with T = 15 K, D = 1 and Ch,i 7 = 1, we there¬ 
fore expect that roughly half of the CO will be depleted once 
the gas density reaches n ~ 15000cm“®. Observationally, the 
density at which xco,gas ~ xco,ice is found to be a factor of 
a few larger than this ( Lippok et al.|[2C)13 l, demonstrating 
that, as expected, our simplified model overestimates the 
effectiveness of CO freeze-out. 

To investigate the effect that freeze-out has on our pre¬ 
dicted CO integrated intensities, we use the simplified model 
outlined above as part of an additional post-processing step. 
Before interpolating the CO densities from our SPH par¬ 
ticles to the Cartesian grid that we use for our radiative 
transfer calculations, we first rescale them by a factor 

77therm 4“ Rcr 


Eaas — 


^CO,gas 

Xco,gas 4- XCO,ii 


77therm 4“ Rcr + R a 


( 12 ) 


to account for the effects of freeze-out, using the expressions 
for Race, 77therm and Rcr given in Equations [7}|^ above. We 
then proceed as before, computing the CO emission from 
the cloud using radmc-3d, and then using the mean values 
of the CO integrated intensity and H 2 column density to 
compute Xco- 

We have carried out this analysis for three of our cloud 
models, Zl-Gl, Z02-G1 and Z003-G1, chosen to span the 
full range of metallicities examined in this study. We list the 
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Table 5. Comparison of Xco computed with and without our 
approximate freeze-out model 


ID 

Xco 



Without freeze-out 

With freeze-out 

Zl-Gl 

1.33 

1.36 

Z02-G1 

3.13 

3.32 

Z003-G1 

21.9 

28.5 


Xco is given in units of the Galactic CO-to-H 2 conversion factor 
^CO,MW (see Eq. 0 and is computed for t = tsp. 

resulting values of Xco in Table aloug with the values 
that we obtain in the absence of freeze-out. We see that 
at solar metallicity, freeze-out has a negligible effect on the 
predicted value of Xco- This is easy to understand, as by the 
time we reach the densities at which CO freeze-out becomes 
signihcant, the gas is highly optically thick in the 1 —> 

0 line. Because the line is optically thick, removing even a 
large fraction of the CO has little effect on the total emission 
(provided that r remains greater than one), and hence little 
effect on Xco- 

As we reduce the metallicity, we find that freeze-out 
starts to have a greater effect on Xco- Again, this is easy 
to understand: at lower Z, we have much less CO and hence 
the optical depth of the line is much smaller, making it more 
sensitive to the actual CO abundance. On the other hand, 
reducing the dust-to-gas ratio, which is assumed in these 
models to scale with Z, reduces the freeze-out rate, since 
there is less dust for the CO molecules to freeze out onto. 
Consequently, in the lower metallicity models, the density 
at which xco,gas ~ xco,ice is much higher than in run Zl- 
Gl. As a result, freeze-out remains relatively unimportant 
in these models despite the lower optical depths, increasing 
Xco by at most 30%. 

Finally, what about the other simulations studied in 
this paper? In the other runs with Go = 1, we would ex¬ 
pect to recover results intermediate between the Zl-Gl and 
Z003-G1 results, with the precise behaviour depending on 
the metallicity. In the runs with reduced dust-to-gas ratio, 
the effects of freeze-out will be smaller, since 7?acc oc T>. 
In the runs with Go = 10 or 100, we would again expect 
the effects of freeze-out to be smaller, since the dust will be 
warmer, making 7?therm larger. Finally, in the runs with ele¬ 
vated cosmic ray ionization rates, Rcr will be larger, and so 
once again we will have less freeze-out. Therefore, it seems 
safe to conclude that the freeze-out of GO onto dust will 
not strongly affect the results we recover for the behaviour 
of Xco in any of our cloud models. 


3.6.2 Magnetic fields 


Real molecular clouds are known to be magnetized (see 


e.g. Crutcher 20121, but we carry out our simulations us¬ 
ing a purely hydrodynamical approach, in view of the tech¬ 
nical difficulties involved in modeling magnetized gas with 
smoothed particle hydrodynamics (see the lengthy discus- 
Price|[2M2 1. The question therefore arises as to how 


much this is likely to affect our results. 

We know from previous simulations of turbulent molec¬ 


ular clouds (e.g. Padoan & Nordlund 2011 Molina et al. 
20121 that prior to the onset of runaway gravitational col¬ 


lapse, the width of the volume density and column density 
PDFs are sensitive to the strength of the magnetic field. 
Increasing the field strength provides additional resistance 
against compression or rarefaction and leads to a narrower 
volume density PDF. Projection of this into two dimensions 
then leads in turn to a narrower column density PDF. The 
extent to which this will affect Xci and Xco will then de¬ 
pend on the size of the characteristic column density (or 
visual extinction) at which we find the transition from C"*" 
to C or from C to CO. If this column density is comparable 
to or smaller than the mean column density of the cloud, 
then a significant fraction of the gas will be [G l] or CO rich, 
regardless of whether we consider a magnetized or unmagne¬ 
tized cloud. In this case, we would expect the values of Xci 
and Xco to be insensitive to the magnetic field strength. 
Similarly, if the required column density is so large that 
it is reached only in regions that go into runaway gravita¬ 
tional collapse (as in e.g. our low metallicity runs), then the 
amount of C and CO that form in the cloud will be insensi¬ 
tive to the width of the turbulence-dominated density PDF, 
implying that Xci and Xco will also not vary much. The 
only case in which the difference in the width of the PDF 
is likely to have a large effect is if the characteristic col¬ 
umn density is such that it is 1-2 standard deviations above 
the mean in the unmagnetized cloud, but multiple standard 
deviations above the mean in the magnetized cloud, owing 
to the change in the size of the standard deviation. In this 
case, we would go from predicting that the cloud is [C l] and 
CO-rich in the unmagnetized case to predicting that it is 
[Cl] and CO-poor in the magnetized case. However, such a 
situation will rarely be encountered in practice, as magnetic 
fields of the level observed in most GMCs do not make a 
large difference to the width of the column density PDF. 
We therefore expect that accounting for this effect may lead 
to minor changes in our results at intermediate metallicities, 
but that it will not greatly change our main conclusions re¬ 
garding the utility of [C l] and CO as tracers of molecular 
gas. 


The other potentially important effect of the magnetic 
held will be to suppress the gravitational collapse of dense 
regions, if they are magnetically sub-critical. This can have a 
large effect on Xco in low metallicity clouds: if gravitational 
collapse is suppressed, the gas will never become CO-bright, 
and so Xco will remain large. Clouds that are so magnet¬ 
ically dominated that they remain starless for long periods 
will therefore also remain invisible in CO for long periods. In 
this circumstance, [C l] will clearly be a much better tracer 
of the clouds. However, it should also be noted that the im¬ 
portance of this effect is unclear. Measurements of the mag¬ 
netic held strength in dense prestellar cores in local molec¬ 
ular clouds suggest that these cores are generally magneti¬ 
cally supercritical, and hence that gravitational collapse is 


not strongly suppressed | Troland & Crutcher 2008 Crutcher 
|2012p . Whether the same holds true for dense clouds in much 
lower metallicity environments is unknown. 


3.6.3 Star formation and feedback 

Our simulations do not account for the effects of stellar feed¬ 
back. Because of this, we cannot explore how Xci and Xco 
evolve over the entire lifetime of one of our simulated clouds, 
as without feedback, there is nothing to stop the whole of 
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the gas eventually being converted to stars. We have there¬ 
fore focussed primarily on their evolution prior to the onset 
of star formation (i.e. at t < tsp), when stellar feedback is, 
by definition, completely absent. 

However, in Section |3.2| we did investigate how Xco 
and Xci evolved in the interval between t = tsF and 
t = tsF -1-0.5 Myr in a few of our simulations, and it is 
therefore important to note that our results for this regime 
may be affected by our neglect of feedback. The importance 
of feedback during this period will depend on the star forma¬ 
tion rate of the clouds and how much of the dense, CO-bright 
gas is actively forming stars. For example, in run Zl-Gl, the 
star formation efficiency of the cloud at t = tsF + 0.5 Myr 
is only around 0.3%, and much of the CO is located in 
gas which is too diffuse to be actively star-forming. In this 
cloud, we would not expect the evolution of Xco and Xci 
at t > tsF to have yet been strongly affected by stellar feed¬ 
back. On the other hand, in run Z003-G100, the star forma¬ 
tion efficiency of the cloud at t = tsF + 0.5 Myr is around 
7.8%, more than an order of magnitude larger than in the 
solar metallicity case. This difference in star formation effi¬ 
ciency is a consequence of the different thermal structure of 
the two clouds. The gas in the low metallicity, strongly irra¬ 
diated cloud is significantly warmer than the gas in the solar 
metallicity cloud. This means that it takes longer to form a 
dense core which is sufficiently massive to be self-gravitating 
(see Tablej^, but also that once such a core has formed, it is 
able to form a lot of stars relatively rapidly. Because of the 
high instantaneous star formation rate in this cloud, and the 
fact that all of the CO is located in a single dense core, it is 
plausible that if we were to include the effects of stellar feed¬ 
back, much of the CO would be destroyed. In this case, the 
period during which the cloud is CO-bright would be even 
more short-lived that we find in our simulations. However, 
confirmation of this requires simulations that do account for 
the effects of stellar feedback, which are out of the scope of 
our present study. 


4 CONCLUSIONS 

The most important point to take away from this study is 
that the usefulness of [C i] emission as a tracer of H 2 mass in 
low metallicity clouds depends strongly on the evolutionary 
state of the clouds. At early times, prior to the onset of 
gravitational collapse in the cloud, [C l] is a far better tracer 
of H 2 than CO. The [C l]-to-H 2 conversion factor, Xci, is not 
constant - it varies with the metallicity roughly as Xci oc 
Z“^ ®-Z~^ for metallicities in the range 0.1 < Z/Zq < 1.0 - 
but it is far less sensitive to changes in the metallicity than 
CO, which is an extremely ineffective tracer of molecular 
mass in starless clouds at low metallicity. 

Shortly before the onset of star formation, however, CO 
becomes a much better tracer of the molecular mass. At 
t = tsF, both Xco and Aci scale approximately with metal¬ 
licity as X oc Z“^, but Xco is consistently smaller than Xci, 
by a factor ranging from a few to an order of magnitude. 
Since it is much easier to observe the '^^CO J = 1 —>■ 0 line 
in terrestrial conditions than the [C l] 609 fj,m line, the con¬ 
clusion we draw from this is that at redshifts close to zero, 
CO should be preferred to [C l] as a tracer of the H 2 mass 


of star-forming clouds in low metallicity systemsj^At higher 
redshifts, however, the [Cl] 609/rm line becomes much eas¬ 
ier to observe from the ground, and in this case [C l] can 
become competitive with CO as a tracer of molecular gas in 


metal-poor star-forming systems | Papadopoulos, Thi & Viti 


20041. 


The difference in the behaviour of starless and star¬ 
forming clouds can be understood as a consequence of the 
physical conditions required to form CO in low metallicity 
clouds. Substantial amounts of CO form only in gas with 
Av ~ 1 or above. In low metallicity clouds, this gas is found 
only in dense, gravitationally collapsing cores. At early times 
in the evolution of the cloud, its structure is dominated by 
turbulence rather than gravity, and these dense cores are not 
present. Consequently, the mean CO brightness of the cloud 
at this time is extremely small. At later times, however, 
dense cores form due to gravitational instability, leading to 
a dramatic increase in the CO brightness of the cloud. On 
the other hand, the [C l] abundance is far less sensitive to 
the visual extinction, and so the mean [C l] brightness of the 
cloud varies much less as the cloud evolves. 

We have also examined how our results change if we 
relax our assumption that the dus-to-gas ratio I? scales lin¬ 
early with the gas-phase metallicity, and instead allow it to 
drop more rapidly with decreasing Z, as suggested by the 
available observational data. We find that in practice, this 
strengthens our results regarding the relative usefulness of 
CO and [Cl] as molecular gas tracers. Signihcantly decreas¬ 
ing leads to a substantial decrease in the CO content of 
the clouds at early times, and hence a substantial increase 
in Xco- On other other hand, it has little influence on the 
value of Xco that we recover once star formation begins. At 
the same time, reducing V also leads to a decrease in the 
amount of H 2 formed in the cloud, which acts to decrease 
Xco, although this effect only seems to be important if the 
dust-to-gas ratio is very low. 

Finally, we have examined the consequences of allowing 
the cosmic ray ionization rate to increase at the same time 
as the strength of the ISRF, as we would expect to occur if 


both directly track the star formation rate (see e.g. Clark & 
|Glover|2015| ). We find that qualitatively, we recover similar 
results, but some of the quantitative details change. Notably, 
increasing the cosmic ray ionization rate by a factor of 10 or 
100 has only a minor effect on the values of Xco we recover 
for our star-forming clouds (contrary to the prediction of 
Bisbas, Papadopoulos, fc Viti|20l5 1, but can lead to a factor 


of two to three decrease in Xci, improving our prospects for 
using [C l] emission to trace H 2 in these systems. 
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